Metabolic rewiring in MYC-driven medulloblastoma by BET-bromodomain inhibition

Medulloblastoma (MB) is the most common malignant brain tumour in children. High-risk MB patients harbouring MYC amplification or overexpression exhibit a very poor prognosis. Aberrant activation of MYC markedly reprograms cell metabolism to sustain tumorigenesis, yet how metabolism is dysregulated in MYC-driven MB is not well understood. Growing evidence unveiled the potential of BET-bromodomain inhibitors (BETis) as next generation agents for treating MYC-driven MB, but whether and how BETis may affect tumour cell metabolism to exert their anticancer activities remains unknown. In this study, we explore the metabolic features characterising MYC-driven MB and examine how these are altered by BET-bromodomain inhibition. To this end, we employed an NMR-based metabolomics approach applied to the MYC-driven MB D283 and D458 cell lines before and after the treatment with the BETi OTX-015. We found that OTX-015 triggers a metabolic shift in both cell lines resulting in increased levels of myo-inositol, glycerophosphocholine, UDP-N-acetylglucosamine, glycine, serine, pantothenate and phosphocholine. Moreover, we show that OTX-015 alters ascorbate and aldarate metabolism, inositol phosphate metabolism, phosphatidylinositol signalling system, glycerophospholipid metabolism, ether lipid metabolism, aminoacyl-tRNA biosynthesis, and glycine, serine and threonine metabolism pathways in both cell lines. These insights provide a metabolic characterisation of MYC-driven childhood MB cell lines, which could pave the way for the discovery of novel druggable pathways. Importantly, these findings will also contribute to understand the downstream effects of BETis on MYC-driven MB, potentially aiding the development of new therapeutic strategies to combat medulloblastoma.

Medulloblastoma (MB) is the most common malignant tumour of the central nervous system (CNS) in children, arising from embryonal lesions generated by diverse progenitor cell populations during early brain development 1,2 . It was first described in 1926 as one subset of glioma 3 . Since then, advances in molecular genetics ameliorated our understanding of MB, culminating in a consensus definition of four distinct groups: wingless (WNT), sonic hedgehog (SHH), Group 3 and Group 4, which differ in both molecular and clinical characteristics 4,5 . Unlike the SHH and WNT groups, no common signalling pathway driving the disease has been identified in Group 3 and 4 tumours. The latter ones exhibit the worst prognosis of the four subgroups and recurrent MYC amplification or overexpression has been identified as one of the major biomarkers for this highrisk MB patient group with poor clinical outcome 4 . The modest understanding of the molecular mechanisms underpinning tumorigenesis in Group 3 and Group 4 MBs has limited the development of targeted therapies. For moderate and high risk-patients trials includes novel chemotherapeutics (pemetrexed and gemcitabine) after standard chemotherapy and risk-adapted radiotherapy. For recurrent and refractory MB patients two trials combine the use of chemotherapeutics with targeted agents such as prexasertib (NCT04023669), a targeted CHK1/2 inhibitor, and ribociclib (NCT01878617), a cyclin-dependent kinase inhibitor 6 . Current treatments include tumour resection, radiotherapy and chemotherapy. Despite this aggressive multimodal therapy, approximatively 30% of the patients die of disease, while survivors suffer long-term side effects due to the harsh treatments on the developing brain 7 www.nature.com/scientificreports/ To date, several preclinical and clinical studies indicated that BET-bromodomain (BRD) inhibitors (BETis) might be considered as next generation agents for treating MYC-driven MB 4,6 . These are small molecules that specifically inhibit proteins in the BRD family which contain two bromodomains capable of recognising acetylated lysine residues in histone tails and recruiting transcriptional factors to promote targeted gene transcription 9 . Their inhibitors prevent the interaction between the bromodomain and the acetyl group, causing the downregulation of certain genes, including c-MYC 10 . JQ1 was the first BETi to be developed. This is a cell permeable small molecule that promotes differentiation and halts proliferation in cancer cell lines and in diverse murine tumour models [10][11][12][13][14] . Despite its efficacy, JQ1 has short-lasting effects due to its half-life of around 1 h, greatly limiting the possibilities to translate preclinical findings into clinical benefit 11,14,15 . OTX-015 (MK-8628) is a BETi developed in the past decade that targets BRD 2/3/4 16 . OTX-015 has been synthesised from ( +)-JQ1 by substituting the tert-butyl ester with a para-hydroxyacetamide. This structural change resulted in an improved pharmacokinetic profile 17 (with an half-life of 6 h) 16,18,19 which led to its use in several clinical trials including glioblastoma multiforme (GBM), NUT midline carcinoma, triple-negative breast cancer (TNBC), castration-resistant prostate cancer, pancreatic ductal carcinoma and haematological malignancies 20 . OTX-015 causes an early, strong and long lasting reduction of MYC, and active doses were not toxic in vivo 16 . Furthermore, in GBM xenograft models this compound has been shown to pass through the blood brain barrier and preferentially bind to cancerous tissue, providing a strong pharmacologic basis for OTX-015 use in brain tumour therapy 21 . Moreover, OTX-015 prodrugs with significantly increased antitumor activity and reduced toxicity have recently been developed 22 .
The c-MYC (MYC) transcription factor belongs to the MYC-proto onco-protein family, which also includes MYCN and MYCL 23 . Heterodimerization with MAX is an essential prerequisite for MYC-driven transcription and cellular transformation. The MYC/MAX dimer binds Enhancer box (E-box) sequences regulating the expression levels of genes involved in growth regulatory networks and oncogenic signalling pathways, including energy production and metabolism 23,24 . Specifically, MYC activates genes modulating metabolic processes in order to balance nutrient supply and demand and drive cell fate decisions 25 . This plays a crucial role in tumour development, as cancer cells reprogram their metabolism to fulfil their increased bioenergetic needs, to maintain a high rate of macromolecular biosynthesis, essential for cell growth/division, and to prevent toxic Reactive Oxygen Species (ROS) levels 26 .
To date, despite a plethora of evidence reporting that MYC tunes metabolic networks to foster oncogenesis in diverse tumour types 27 , data remain scarce and limited, and likewise, the metabolic response of BETi treatment of MYC-driven MB remains unexplored. In this work, we investigate changes in metabolism of MYC-driven MB cell lines prior and after treatment with OTX-015. We chose OTX-015 over other inhibitors, e.g. JQ1, because of its higher promise in clinical trials. This study contributes to the phenotypical characterisation of MYC-driven MB models and reveals metabolic changes caused by OTX-015, providing new clues of the anticancer mechanisms of action for these compounds.

OTX-015 inhibits cell viability and reduces MYC expression in medulloblastoma cell lines.
We selected human D283 and D458 MB cell lines as in vitro models of MYC-driven MB. The D283 cell line has been classified as either Group 3 or Group 4 and overexpresses MYC and the transcription factor orthodenticle homeobox 2 (OTX2) 28 . The latter plays a key role in the pathogenesis of anaplastic MB and treatment with all trans-retinoic acid causes OTX2 downregulation and apoptosis 29 . D458 cells belong to Group 3 and are characterised by MYC amplification. They derive from the metastasis of a parental tumour and harbour wild type p53 28 with R72P polymorphism, a common single nucleotide polymorphism (SNP) that results in either an arginine or proline at position 72 of the protein 30,31 .
The D283 and D458 cells were treated with a range of OTX-015 doses spanning from 0 to 10 µM for 24 and 48 h, to evaluate the optimal time point to study their metabolism. Western Blot analysis was performed to evaluate MYC protein expression ( Fig. 1A-C and Supplementary Fig. S1). Using the MTT assay we demonstrated that OTX-015 decreases cell viability by 50% after 48 h, while Western Blot showed that c-MYC protein levels were already reduced after 24 h of incubation with OTX-015. Specifically, relative protein quantification showed that 7.5 µM of OTX-015 caused a 70% decrease in c-MYC protein levels versus untreated cells. We thus choose to investigate the early effects of BET inhibition and MYC downregulation on cancer metabolism after 24 h of incubation with 7.5 µM of OTX-015.
OTX-015 treatment causes a statistically significant metabolic shift in MYC-driven MB cell lines. Next, we analysed the metabolism of D283 and D458 cells upon OTX-015 treatment by Nuclear Magnetic Resonance (NMR) spectroscopy. 1 H NMR profiling of the aqueous cellular extracts allowed the identification of 41 intracellular polar metabolites (Fig. 2, Supplementary Figs. S2-S7, and Supplementary Table S1), assigned by comparing peak chemical shifts, J-couplings and multiplet intensity ratios to NMR spectral resonance patterns available in the Human Metabolome Database (HMDB) 32 and the biological magnetic resonance data bank (BMRB) 33  www.nature.com/scientificreports/ amino acid derivatives including reduced glutathione GSH, creatine, taurine, hypotaurine and formate; phospholipid derivatives including phosphocholine and glycerophosphocholine; nucleotides as ATP, UDP, CTP and NAD; the carbohydrate myo-inositol and its stereoisomer scyllo-inositol; the polyamine putresceine and the vitamin B5 pantothenate. Furthermore, two phosphorylated sugars were observed to change upon OTX-015 treatment. One of these was identified based on the chemical shift of the anomeric proton as UDP-N-acetylglucosamine, whilst the other was dubbed as a generic 'UDP-sugar' , as the unambiguous assignment of its structure was not possible from our NMR data alone. 1D 1 H NMR spectra were used to perform a Partial least squares discriminant analysis (PLS-DA), a supervised multivariate statistical analysis which is a powerful tool for dealing with large-scale metabolomics studies 34 . The PLS-DA score scatter plot shows that the metabolome of D283 and D458 cell lines is statistically different, with component 1 and component 2, with a variance of 46.4% and 12% respectively, responsible for group separation (Fig. 3A). Notably, a clear change of the metabolome is detectable in both cell lines upon the treatment with OTX-015, as D283_OTX and D458_OTX are shifted along component 2 compared to their respective controls (D283_Ctrl and D458_Ctrl). Next, to identify the metabolites that contribute to the observed PLS-DA separation, variable importance in projection (VIP) scores were determined for the NMR data set (Fig. 3B). VIP scores ˃1 suggest that dissimilar levels of glycine, taurine, myo-inositol, glutamate, creatine, and N-acetyl aspartate are strongly associated with the separation of D283 and D458 cells, whereas increased myo-inositol and taurine, and reduced glycine, glutamate and N-acetyl aspartate determine PLS-DA separation between treatment and controls.

Key metabolites regulated by OTX-015 are shared in both MB cell lines. A heat map representing
a quantitative targeted analysis confirmed key metabolic differences between the groups identified by PLS-DA (Fig. 4A). One-way analysis of variance (ANOVA) for this metabolite quantification is presented in Table 1. Fisher's Least Significant Difference (LSD) indicates that D283 and D458 are characterised by different levels of amino acids and derivatives (N-acetylaspartate, taurine, hypotaurine, tyrosine, glycine, phenylalanine, serine, leucine, isoleucine, valine, alanine, GSH, glutamate, and threonine); carbohydrates (myo-inositol, scyllo-inositol, UDP-N-acetylglucosamine, and UDP-sugar); phospholipid derivatives (phosphocholine and glycerophosphocholine); nucleotides (NAD and CTP); the polyamine putresceine, and the vitamin B5 pantothenate.
Upon treatment with OTX-015, a metabolic response was detected in both D283 and D458 cell lines and interestingly a few metabolites appear to be commonly regulated. These include myo-inositol, glycerophosphocholine, serine, pantothenate, UDP-N-acetylglucosamine, phosphocholine and glycine, all of which increased in both D283 and D458 after treatment ( Fig. 4B-D). Notably, OTX-015 treatment induced additional changes in the endogenous metabolome of D283 and D458, albeit the specific alterations differ between the two cell lines ( Fig. 4A and Supplementary Fig. S11).  Table 1), along with their relative concentrations, using the web-based platform MetaboAnalyst 35 (Fig. 5). The top 25 pathways significantly enriched in D283_OTX versus untreated control and D458_OTX versus untreated control cells are shown in Fig. 5A,B. We focused on the common enriched pathways in which the metabolites myo-inositol, glycerophosphocholine and serine increased significantly upon OTX-015 treatment in both cell lines (Fig. 5). These data confirm that OTX-015 triggers a different metabolic reprogramming in D283 and D458 MB cell lines, reflecting their diverse basal metabolism. Nonetheless, in both cell lines OTX-015 affects ascorbate and aldarate metabolism, inositol phosphate metabolism, phosphatidylinositol signalling system, ether lipid metabolism, glycerophospholipid metabolism, glycine serine and threonine metabolism and aminoacyl tRNA biosynthesis (Fig. 5C).

Discussion
MYC-driven MB is an aggressive prevalent paediatric tumour, which is refractory to intensive multimodal therapy 4,36 . In this scenario, BETis are emerging as promising agents in MYC-dependent MB as well as in other tumours, because of their ability to impair MYC expression and to modulate the epigenetic machinery of cancer cells 4,20 . In this work, we studied the effect of the BETi OTX-015 on the metabolism of the D283 and D458 MYCdriven childhood MB cell lines. OTX-015 has been developed from the conventional BETi JQ1 to improve its pharmacokinetic profile and thus its potential in clinical use 16,17,22 .
It is already known that MYC regulates metabolism which in turn controls MYC levels 27 . MYC promotes cell proliferation and induces metabolic alterations because of increased anabolic demand. On the other hand, MYC is also downstream of multiple control mechanisms regulated by nutrient levels and metabolic stress 27 . Despite the well-recognised crosstalk between MYC and metabolism in cancer, how metabolism is deregulated in MYC-driven MB is still unknown. Furthermore, although inhibition of BRD proteins is known to cause metabolic reprogramming 37,38 , whether BETi-dependent metabolic changes can be used to monitor drug response towards BETis is not understood. Here, we investigated the metabolome of MYC-driven MB cell models and interrogated how OTX-015 would modulate their metabolic phenotype using an NMR-based approach. NMR spectroscopy has proven to be of high value as an accurate and quantitative method to investigate metabolic alterations related to disease, through detecting metabolic signatures for early diagnosis and prognosis evaluation, unveiling cellular mechanisms and discovering pharmacological actions and potential therapeutic significance of drug compounds [39][40][41] .
In this work, we discovered that D283 exhibits diverse metabolic characteristics compared to D458 cells, consistent with their different origin and genetic background 28 . The most significantly different metabolites between the two cell lines were: amino acids and their derivatives such as hypotaurine and GSH; the TCA cycle metabolite citric acid; the phospholipid derivatives phosphocholine and glycerol-phosphocholine; the vitamin B5 pantothenate; NAD and scylloinositol. Upon treatment of D283 and D458 with OTX-015, common features were revealed in addition to cell-specific metabolic changes. In both cell lines, quantitative NMR data showed www.nature.com/scientificreports/ an increase in myo-inositol, a metabolite well known to control cell osmoregulation, whose role in pathophysiological contexts is tissue-and tumour-dependent 42 . Notably, myo-inositol targets critical cancer liabilities [43][44][45] by: (i) inhibiting pRB phosphorylation and subsequently blocking the cell cycle; (ii) reducing PI3K, Akt, ERK and NF-kB and (iii) impairing invasiveness [43][44][45] . The observed increase in myo-inositol levels therefore may contribute to the anti-proliferative effect observed in D283 and D458 cell lines upon treatment with OTX-015. Metabolite set enrichment analysis (MSEA) revealed that OTX-015 treatment causes a marked common enrichment in metabolic pathways involving myo-inositol (enrichment ratio ˃8). These include ascorbate and aldarate metabolism, inositol phosphate metabolism and phosphatidylinositol signalling system pathways in both cell lines. To our knowledge, ascorbate and aldarate metabolism has never previously been linked to MB, albeit it has been described as markedly altered in other cancers 46 . Notably, both ascorbate (AA) and its oxidised form dehydroascorbate (DHA), which together represent the two existing redox states of Vitamin C, are important   48 . In agreement with previous work 49 , our results suggest that, beyond MYC downregulation, OTX-015, and perhaps other BETis, could elicit their antitumour capabilities also by disrupting the metabolic equilibrium between antioxidants and ROS, thereby eventually leading to cancer cell death.
Our data also show that inositol phosphate metabolism and phosphatidylinositol signalling system pathways are associated with OTX-015 treatment. These two pathways are well known to be initiated upon stimulation www.nature.com/scientificreports/ of phospholipase C, which metabolises phosphatidylinositol 4,5-bisphosphate into the intracellular second messengers 1,2-diacylglycerol (DAG) and inositol 1,4,5-trisphosphate (IP3) 50 . DAG activates protein kinase C, whilst IP3 regulates cellular calcium, and both can be converted into more polar inositol phosphates (IPs) and diphosphoryl inositol phosphates (PP-IPs). Inositols, of which the most abundant stereoisomer is myo-inositol 51 , play a crucial role in several biological functions and deregulation of their synthesis has been linked to disease including cancer 43,44 . IPs and PP-IPs are key in transcription control, mRNA export and DNA repair, and are also important for glucose sensing in β-pancreatic cells 52 . As mentioned above, inositols can impair oncogenic signalling using several mechanisms including inhibition of the PI3K/Akt pathway 43 . Notably, this has already been recognised in medulloblastoma as one of the relevant oncogenic pathways contributing to the disease 53 and, accordingly, a selective pan-class I PI3K inhibitor namely BKM120 has been revealed as a potential novel therapy for MB 54 . Thus, in our MYC-driven cell models OTX-015 treatment may upregulate myo-inositol production, counteracting PI3K/Akt/mTOR signalling and eventually triggering tumour cell death. www.nature.com/scientificreports/ Our NMR analysis recorded higher concentrations of glycerophosphocholine and phosphocholine upon OTX-015 treatment. Increase in these metabolites has been strongly associated with the well-known accumulation of lipid droplets 55,56 . Previous work demonstrated that MYC inhibition, achieved either with the c-MYC/ Max inhibitor 10058-F4 or by MYC depletion, triggers mitochondrial dysfunction and β-fatty acid oxidation impairment, eventually leading to accumulation of intracellular lipid droplets 57,58 . Notably, BRD4 inhibition has also been shown to regulate OXPHOS capacity and to rewire metabolism in complex I (CI) mutant cybrid cells 37 . Consistently with this, our data appear to suggest that one of the putative anticancer mechanisms of OTX-015 may be hindering mitochondrial function, which would lead to inefficient lipid catabolism and increased levels of glycerophosphocholine and phosphocholine. Further experiments are however needed to investigate the intriguing possibility of a crosstalk between BETis, the mitochondrial machinery and lipid metabolism in MYC-driven MB.
In agreement with the increased levels of glycerophosphocholine, our MSEA analysis revealed that OTX-015 also causes a common enrichment of two intertwined metabolic pathways involving this metabolite, specifically ether lipid metabolism and glycerophospholipid metabolism. Although diverse roles of these pathways in cancer have been reported [59][60][61][62] , certain types of ether lipids function as endogenous antioxidants and are involved in cell differentiation 63 , whilst glycerophospholipids have been found decreased in glioma stem-like cells 64 . Changes in these two pathways may therefore have an impact on the grade of cell differentiation of MYC-driven MB.
Upon treatment with OTX-015, we also detected a significant concentration rise of serine in both cell lines. This amino acid is synthesised from the glycolytic intermediate 3-phosphoglycerate and can be transformed into glycine, which in turn is used in one carbon metabolism processes to sustain the synthesis of macromolecule building blocks such as purines, and/or for GSH synthesis to counteract ROS effects 65 . Importantly, serineglycine conversion is catalysed by the serine hydroxymethyltransferase (SHMT), of which both the cytoplasmic (SHMT1) and mitochondrial (SHMT2) isoforms are known to be upregulated by MYC 66 . Thus, we could speculate that downregulation of MYC caused by OTX-015 impairs SHMT1/2 activity, resulting in serine accumulation. Consistent with this, we also detected a concomitant reduction in glycine in both cell lines and revealed that a commonly affected pathway following OTX-015 treatment is glycine, serine and threonine metabolism. Importantly, glycine is a precursor of purines and thus a key metabolite for rapidly proliferating cancer cells 67 . Accordingly, a recent study unveiled that glycine levels are increased in MYC-amplified MB with respect to the normal cerebellum 68 .
Another pathway that is affected in both cell lines in response to OTX-015 treatment is aminoacyltRNA biosynthesis. Notably, supporting our results, aminoacyltRNA biosynthesis has been reported as a pathway signature for group 3γ MB tumours, which is characterised by increased MYC copy number and negative prognosis 69 . AminoacyltRNA synthases are key enzymes in protein synthesis that during evolution acquired further domains, which provided these enzymes with additional roles beyond translation. Upregulation of seryl-tRNA synthetase induced cellular senescence and inhibited the growth of cervical tumour xenografts in mice by causing the senescence of tumour cells 70,71 . We found that MYC-driven MB cells treated with OTX-015 have significant higher levels of serine compared to untreated cells, and this is compatible with an upregulation of seryl-tRNA synthetase 70 . OTX-015 may therefore promote senescence by enhancing seryl-tRNA synthetase biosynthesis, although any link between OTX-015 and seryl-tRNA biosynthesis awaits to be proven.
In summary, in this work we describe the metabolic profile of MYC-driven MB cell models and reveal key metabolic changes triggered by OTX-015 treatment that could be either MYC-dependent or -independent. These novel insights contribute to the phenotypical characterisation of MYC-driven MB and could lay the foundation for uncovering new druggable signalling pathways beyond MYC-inhibition. Notably, our analysis revealed how OTX-015 was able to alter the metabolism of MYC-driven MB cell models. This provides key additional information on possible mechanisms of actions elicited by this compound at the cellular level to inhibit the growth of susceptible cancer cells, whilst also generating potentially useful clues for the use of BETis in the treatment of medulloblastoma. This first metabolic characterisation is likely to promote, guide and pave the ways to further analyses, which will likely reveal additional metabolic features and possible ramifications for the medulloblastoma treatment with OTX-015 or other BETis.
Western blot. Whole-cell lysates were prepared using RIPA buffer (Sigma-Aldrich) supplemented with Halt Protease and Phosphatase Inhibitor Mixture (Thermo Scientific). Samples were resolved on SDS/PAGE and transferred to Polyvinylidene fluoride (PVDF) membranes (BioRad, 1704156). The c-Myc antibody (9E10) (NB600-302, Novus Biologicals Bio-Techne, United Kingdom 1:2000) was used to identify c-MYC, while anti-βactin (Santa Cruz, SC47778 1:10,000) was used as loading control. HRP-conjugated secondary antibodies were from GE Healthcare. Bands were detected using SuperSignal West Dura chemiluminescent substrate (Thermo Fisher Scientific, 34075) and images were obtained in a ChemiDoc XRS + (BioRad). Image processing was performed with ImageLab (Life Science Research, BioRad). All Western blot data are the result of three independent biological replicates. The quantification was performed by normalising to β-actin and the percentage of decrease in c-MYC in treated versus control cells was calculated.  72 , which was then phase corrected, baseline corrected. The chemical shifts were referenced to the TSP peak at 0.0 ppm. Spectra were then exported to MATLAB, using a toolbox written by A. Edison and co-workers 73 . Spectra were aligned by the Pearson PAFFT method (peak alignment by fast Fourier transform) and then normalised by probabilistic quotient normalisation (PQN).
Six biological replicates were used for this analysis.
Metabolite identification and quantification. 2D spectra were used on one of the samples for identification purposes, consisting of a Total Correlation SpectroscopY ([ 1 H-1 H] TOCSY), using the Bruker pulse sequence "dipsi2gpphzs", slightly modified to include presaturation, with eight scans, 256 t 1 increments, a spectral width of 13.7 ppm in both dimensions and a relaxation delay of 2 s, along with a Heteronuclear Single-Quantum Correlation (([ 1 H-13 C] HSQC) spectrum, using the Bruker pulse sequence "hsqcetgpsisp", with 16 scans, 256 t 1 increments, a spectral width of 170 ppm in the 13 C dimension and 12 ppm in the 1 H dimension, and a relaxation time of 2 s. Metabolites were assigned using information acquired from the human metabolome database (HMDB) 32 and the biological magnetic resonance data bank (BMRB) 33 . Relative quantification of metabolites was achieved by integration of peaks to calculate peak area, which is proportional to the number of nuclei that give rise to the signal. Mean ± SD values of peak areas were calculated for each metabolite.
Statistical analysis using MetaboAnalyst. After identification of the metabolites, a list of representative peaks was chosen for each of the identified metabolites (Supplementary Table 1). Their peak area for each spectrum was taken and stored in a table which, as a .csv file, was then loaded onto MetaboAnalyst 35 for the statistical analysis. For multivariate analyses, spectra were scaled by the pareto method 74 . To determine whether spectra cluster into groups, PCA and PLS-DA with k-fold cross-validation were performed, along with a one-way ANOVA to see which metabolites were statistically different for specific group comparison.
Metabolic pathways enrichment analysis using MetaboAnalyst. The signal integrals corresponding to the identified metabolites were used to perform a Quantitative Enrichment Analysis employing MetaboAnalyst 5.0 35 platform with the KEGG library, a database of relevant metabolic pathways found in humans [75][76][77] . Significant pathways were those with p-value < 0.05. The top 25 enriched pathways were shown.

Data availability
A large part of data generated or analysed during this study are included in this published article and its supplementary information files. Some of the datasets used and/or analysed during the current study have not been included and are available from the corresponding authors on reasonable request.